#include <assert.h>
#include <stdio.h>

#include "../src/mstep.c"

typedef struct _MstepTestDatum MstepTestDatum ;
struct _MstepTestDatum 
{ 
  struct
    {
      float    loga[6 * 4];
      float    logbeta;
      float    c[6 * 5];
      float    d[6 * 5];
      float    pz[5];
      Rating   ratings[8];
      float    logbetamean;
      float    logbetastddev;
    }                   inputs;

  struct
    {
      float    newq;
      float    logbeta;
    }                   desired_output;
};

static void
test_mstep (void)
{
  MstepTestDatum data[] = 
    {{{{0.953769, 0.861753, 0.201633, 0.791226, 0.128144, 0.504315,    0.346759, 0.0821954, 0.442096, 0.000539102, 0.247704, 0.829286,    0.971198, 0.607911, 0.43943, 0.329488, 0.396738, 0.697844,    0.684552, 0.817889, 0.745096, 0.574393, 0.203272, 0.1696},   0.791328, {0.71264, 0.00163917, 0.378373, 0.663184, 0.208325,    0.65488, 0.296178, 0.221088, 0.207786, 0.407176, 0.466892,    0.249889, 0.599875, 0.967746, 0.137404, 0.853151, 0.902031,    0.283194, 0.319514, 0.108055, 0.327637, 0.0799223, 0.149915,    0.316727, 0.614997, 0.0782831, 0.771542, 0.653543, 0.406673,    0.423403}, {0.475364, 0.432456, 0.198887, 0.0162273, 0.00847145,    0.182567, 0.599012, 0.0484811, 0.871068, 0.329415, 0.696982,    0.765287, 0.551553, 0.22136, 0.369344, 0.685365, 0.401638,    0.904633, 0.754347, 0.607082, 0.630097, 0.25109, 0.347674,    0.183678, 0.154733, 0.818634, 0.148788, 0.167451, 0.146262,    0.636067}, {0.549776, 0.11897, 0.275194, 0.306652,    0.852794}, {{3, 0}, {3, 1}, {2, 2}, {3, 2}, {1, 4}, {2, 2}, {1,     0}, {3, 0}}, 0.353683,   0.723641}, {-30.393941468260973978, -1.2932330518308693577}}, {{{0.0852915, 0.48345, 0.668318, 0.322003, 0.180658, 0.729103,    0.0612366, 0.691906, 0.929569, 0.381428, 0.877558, 0.537173,    0.110935, 0.232641, 0.710107, 0.390912, 0.474868, 0.682865,    0.591138, 0.115718, 0.168216, 0.830071, 0.237455, 0.392077},   0.0829245, {0.346621, 0.569137, 0.0700739, 0.902266, 0.617519,    0.5079, 0.378168, 0.972697, 0.23609, 0.630342, 0.840994, 0.861762,     0.00344937, 0.920234, 0.450083, 0.386894, 0.320584, 0.329097,    0.334365, 0.218678, 0.490513, 0.0916417, 0.942289, 0.135754,    0.143892, 0.522505, 0.872215, 0.233488, 0.526373,    0.0146049}, {0.494047, 0.26079, 0.290283, 0.384263, 0.653053,    0.399028, 0.286834, 0.464029, 0.20297, 0.0121343, 0.966249,    0.134933, 0.868604, 0.793456, 0.475736, 0.0432909, 0.926316,    0.657703, 0.331845, 0.520786, 0.054101, 0.424215, 0.805472,    0.506181, 0.560054, 0.163425, 0.515189, 0.121918, 0.907001,    0.764396}, {0.228355, 0.657889, 0.704032, 0.752262,    0.262106}, {{3, 1}, {1, 1}, {2, 1}, {0, 1}, {1, 4}, {3, 1}, {3,     1}, {0, 4}}, 0.522956,   0.835427}, {-36.454907382172583304, -1.3864568431145718404}}, {{{0.958806, 0.78637, 0.479665, 0.909111, 0.301103, 0.454525,    0.958879, 0.85501, 0.876888, 0.649053, 0.452698, 0.294956,    0.713463, 0.133864, 0.330781, 0.387955, 0.949067, 0.905509,    0.672892, 0.683924, 0.196805, 0.643403, 0.149936, 0.848496},   0.238, {0.857034, 0.670271, 0.939385, 0.936897, 0.402509, 0.711391,     0.0843745, 0.0600087, 0.753456, 0.258693, 0.789418, 0.346545,    0.619591, 0.927912, 0.401463, 0.397478, 0.714082, 0.25502,    0.717539, 0.200672, 0.0706789, 0.105084, 0.869043, 0.962673,    0.213645, 0.434813, 0.929658, 0.025776, 0.811136,    0.723422}, {0.845284, 0.965767, 0.0576805, 0.464729, 0.0558656,    0.619222, 0.438089, 0.536817, 0.654403, 0.221744, 0.724007,    0.281796, 0.936863, 0.0210717, 0.653328, 0.176712, 0.0678203,    0.058399, 0.439683, 0.741899, 0.138162, 0.0326231, 0.628546,    0.0184766, 0.292879, 0.0668557, 0.570866, 0.553748, 0.237013,    0.447634}, {0.132777, 0.0169312, 0.58261, 0.225889,    0.40877}, {{2, 3}, {2, 1}, {1, 1}, {3, 0}, {0, 0}, {3, 1}, {1,     4}, {1, 0}}, 0.735135,   0.645747}, {-24.335248771610502524, -1.2344815572064789905}}, {{{0.204818, 0.755442, 0.558423, 0.577927, 0.146419, 0.315759,    0.816524, 0.439765, 0.113796, 0.687213, 0.798048, 0.146886,    0.0469398, 0.116347, 0.2443, 0.909873, 0.599306, 0.98357,    0.227369, 0.327262, 0.373417, 0.5748, 0.492234, 0.681515},   0.168599, {0.819358, 0.933811, 0.103589, 0.0221805, 0.503599,    0.117287, 0.663824, 0.908385, 0.816386, 0.319239, 0.516938,    0.861445, 0.700039, 0.0749392, 0.607065, 0.262139, 0.716469,    0.84757, 0.279803, 0.888722, 0.141669, 0.355336, 0.598287,    0.720123, 0.322311, 0.421525, 0.494699, 0.697942, 0.818712,    0.304239}, {0.830875, 0.789557, 0.00232642, 0.985, 0.313937,    0.928112, 0.302287, 0.91006, 0.706871, 0.665974, 0.585818,    0.0624901, 0.427069, 0.777252, 0.444149, 0.707154, 0.828781,    0.0571287, 0.121838, 0.285628, 0.334083, 0.359186, 0.303125,    0.981389, 0.503208, 0.569629, 0.300799, 0.99639, 0.189272,    0.641516}, {0.998511, 0.0863294, 0.4824, 0.975543,    0.412693}, {{1, 2}, {0, 3}, {3, 3}, {3, 1}, {1, 2}, {2, 3}, {2,     4}, {2, 2}}, 0.0238393,   0.0553316}, {-63.852194532909647822, -0.12744441219547256339}}, {{{0.198291, 0.968545, 0.316686, 0.22655, 0.141163, 0.846707,    0.0310577, 0.892467, 0.781976, 0.543582, 0.0496682, 0.389259,    0.212348, 0.242783, 0.0532784, 0.199988, 0.570831, 0.244272,    0.966949, 0.717587, 0.595288, 0.831579, 0.94311, 0.662256},   0.396997, {0.863034, 0.626424, 0.435706, 0.255834, 0.0163266,    0.595366, 0.543238, 0.473858, 0.472744, 0.545698, 0.153979,    0.261511, 0.229961, 0.49242, 0.953992, 0.69068, 0.985689,    0.525471, 0.236404, 0.0953912, 0.15411, 0.582361, 0.574148,    0.698394, 0.291077, 0.955937, 0.138443, 0.44256, 0.27475,    0.360571}, {0.595205, 0.968702, 0.802006, 0.814873, 0.441225,    0.707191, 0.572045, 0.322453, 0.487234, 0.0165113, 0.586356,    0.796983, 0.25083, 0.92112, 0.432245, 0.214622, 0.676681,    0.222726, 0.141169, 0.258685, 0.538238, 0.780166, 0.866419,    0.898114, 0.943034, 0.811465, 0.0644134, 0.0832407, 0.501808,    0.104274}, {0.492369, 0.760787, 0.0145745, 0.0877625,    0.906013}, {{1, 3}, {2, 0}, {2, 2}, {1, 2}, {1, 0}, {2, 1}, {3,     0}, {1, 3}}, 0.963805,   0.763745}, {-36.172573614166577558, -1.4048061843050376653}}, {{{0.166642, 0.473768, 0.749183, 0.0870637, 0.943916, 0.332599,    0.490498, 0.548825, 0.16375, 0.46618, 0.592384, 0.605792,    0.352285, 0.401767, 0.509144, 0.103983, 0.248012, 0.909398,    0.748356, 0.0894086, 0.160249, 0.00338457, 0.784552, 0.325664},   0.993607, {0.529617, 0.0353691, 0.2386, 0.0496903, 0.197018,    0.544871, 0.689775, 0.88594, 0.730838, 0.952487, 0.0839833,    0.533655, 0.329071, 0.443343, 0.98, 0.285643, 0.419673, 0.694986,    0.890592, 0.125394, 0.416289, 0.910434, 0.564928, 0.131787,    0.886672, 0.875065, 0.326328, 0.0820969, 0.689654,    0.330194}, {0.636553, 0.196157, 0.958817, 0.377708, 0.55257,    0.662502, 0.629745, 0.934365, 0.572569, 0.376859, 0.210072,    0.239379, 0.681978, 0.251465, 0.793783, 0.328944, 0.11705,    0.119678, 0.907111, 0.453879, 0.790722, 0.0375813, 0.217457,    0.123685, 0.15417, 0.841424, 0.25864, 0.745977, 0.6016,    0.178922}, {0.628895, 0.811612, 0.0290307, 0.802063,    0.418823}, {{1, 2}, {0, 4}, {1, 4}, {0, 2}, {0, 3}, {2, 2}, {3,     4}, {0, 3}}, 0.572233,   0.347053}, {-49.366051075610333736, -0.90186748935132083315}}, {{{0.550598, 0.62504, 0.243289, 0.230003, 0.43092, 0.717929, 0.78941,    0.439281, 0.393339, 0.500472, 0.665726, 0.285111, 0.551914,    0.241832, 0.919749, 0.683511, 0.372992, 0.612938, 0.108137,    0.65448, 0.570929, 0.194115, 0.535903, 0.307427},   0.0203306, {0.569075, 0.292614, 0.0774244, 0.589411, 0.851146,    0.503204, 0.638144, 0.196072, 0.350673, 0.837478, 0.353033,    0.644158, 0.108841, 0.917729, 0.669522, 0.271166, 0.495903,    0.809593, 0.0150413, 0.700238, 0.301789, 0.273689, 0.707614,    0.679907, 0.732714, 0.981075, 0.63019, 0.0904965, 0.881568,    0.477872}, {0.992046, 0.894424, 0.530895, 0.640394, 0.639013,    0.250266, 0.422054, 0.722664, 0.969491, 0.9791, 0.92615, 0.913072,     0.95445, 0.278862, 0.624361, 0.639382, 0.246836, 0.598955,    0.891647, 0.658307, 0.616647, 0.508458, 0.0100793, 0.180435,    0.624601, 0.614034, 0.479185, 0.540042, 0.985588,    0.363768}, {0.057131, 0.817377, 0.0160962, 0.384668,    0.130981}, {{1, 2}, {1, 1}, {2, 2}, {1, 0}, {2, 4}, {2, 3}, {3,     0}, {0, 1}}, 0.904305, 0.0616461}, {-78.716847843441618507,   0.63185047284487011473}}, {{{0.105806, 0.50662, 0.264923, 0.81481,    0.506852, 0.614972, 0.606616, 0.198163, 0.998393, 0.604893,    0.426181, 0.573563, 0.384359, 0.125708, 0.886139, 0.587975,    0.0205915, 0.0685771, 0.068762, 0.571879, 0.635923, 0.937596,    0.164457, 0.510233},   0.530117, {0.430977, 0.899534, 0.695423, 0.0232656, 0.816005,    0.292918, 0.497259, 0.0248723, 0.211112, 0.866738, 0.923697,    0.640513, 0.0854036, 0.980599, 0.335722, 0.619921, 0.0168265,    0.911837, 0.763843, 0.983998, 0.0792303, 0.74738, 0.25361,    0.453881, 0.648254, 0.847846, 0.558188, 0.430615, 0.832249,    0.554928}, {0.0609282, 0.405743, 0.621137, 0.68819, 0.137231,    0.76523, 0.535734, 0.707592, 0.80151, 0.145309, 0.518907,    0.795755, 0.0376666, 0.161311, 0.439677, 0.0483752, 0.784056,    0.70743, 0.791423, 0.200529, 0.225869, 0.276815, 0.959174,    0.645601, 0.16494, 0.871072, 0.338037, 0.957411, 0.0277089,    0.105842}, {0.802303, 0.24982, 0.226199, 0.960533,    0.283396}, {{0, 2}, {1, 3}, {3, 3}, {0, 3}, {3, 4}, {1, 1}, {3,     3}, {3, 2}}, 0.454065,   0.188533}, {-55.146760336174218499, -0.47683249966028373746}}, {{{0.799222, 0.843719, 0.405689, 0.404476, 0.0917926, 0.0522962,    0.20516, 0.178608, 0.814978, 0.093122, 0.559559, 0.0136674,    0.943907, 0.755085, 0.602147, 0.985959, 0.838065, 0.952781,    0.352328, 0.759759, 0.877532, 0.669385, 0.898263, 0.571227},   0.0783098, {0.825666, 0.492574, 0.16675, 0.986517, 0.773369,    0.287414, 0.988143, 0.171539, 0.680247, 0.727855, 0.974475,    0.227633, 0.925163, 0.125708, 0.988517, 0.389568, 0.972381,    0.77338, 0.228757, 0.512035, 0.302996, 0.875117, 0.657531,    0.433726, 0.47733, 0.382543, 0.49078, 0.447208, 0.703961,    0.0951289}, {0.502638, 0.275669, 0.0237133, 0.367274, 0.528163,    0.0480369, 0.0985508, 0.241566, 0.539646, 0.658469, 0.12617,    0.468186, 0.310889, 0.146434, 0.823174, 0.59307, 0.653358,    0.712708, 0.345844, 0.210527, 0.162578, 0.2655, 0.641883,    0.115398, 0.65994, 0.98983, 0.618169, 0.748124, 0.131778,    0.941793}, {0.519619, 0.506558, 0.592132, 0.283324,    0.393449}, {{0, 3}, {0, 0}, {3, 1}, {2, 3}, {1, 0}, {3, 0}, {1,     1}, {0, 3}}, 0.0383717,   0.281243}, {-37.226939711783911656, -0.84433872332147009936}}, {{{0.13689, 0.570275, 0.445302, 0.627884, 0.424182, 0.224432,    0.234775, 0.465306, 0.158682, 0.582549, 0.119377, 0.805366,    0.168852, 0.964379, 0.371253, 0.673588, 0.227058, 0.444761,    0.864695, 0.0814564, 0.943734, 0.0513116, 0.826323, 0.800214},   0.806844, {0.481036, 0.381021, 0.172329, 0.382662, 0.256605,    0.146246, 0.707023, 0.22398, 0.674056, 0.026869, 0.901657,    0.0551281, 0.709677, 0.655616, 0.228069, 0.82807, 0.264916,    0.790922, 0.146613, 0.884335, 0.213604, 0.964599, 0.346399,    0.0774913, 0.732568, 0.583578, 0.17407, 0.694829, 0.475963,    0.437332}, {0.467047, 0.470849, 0.801907, 0.410463, 0.565389,    0.415721, 0.0922308, 0.754847, 0.33732, 0.587651, 0.827315,    0.963925, 0.190707, 0.703316, 0.613711, 0.999326, 0.844308,    0.625825, 0.881143, 0.415748, 0.670238, 0.930996, 0.405179,    0.978416, 0.203191, 0.460146, 0.603272, 0.567953, 0.637802,    0.0444251}, {0.511041, 0.813106, 0.300482, 0.456774,    0.683726}, {{3, 2}, {2, 3}, {3, 3}, {3, 1}, {3, 0}, {1, 2}, {3,     0}, {1, 3}}, 0.849181,   0.109775}, {-105.92714695213845007, -0.0051218544519903359743}}};

  for (unsigned int i = 0; i < sizeof (data) / sizeof (data[0]); ++i)
    {
      NormalDistribution prior_logbeta = 
        normal_distribution (data[i].inputs.logbetamean, 
                             data[i].inputs.logbetastddev);

      float hq = 0;
      for (unsigned int k = 0; k < 5; ++k)
        {
          hq -= data[i].inputs.pz[k] * log (data[i].inputs.pz[k]);
        }

      float logbeta = data[i].inputs.logbeta;
      float q = mstep (data[i].inputs.loga,
                       6,
                       &logbeta,
                       data[i].inputs.c,
                       5,
                       data[i].inputs.d,
                       data[i].inputs.pz,
                       data[i].inputs.ratings,
                       8,
                       &prior_logbeta.base,
                       hq,
                       1000,
                       1e-7);

      assert (fabs (q - data[i].desired_output.newq) <=
              1e-4 * (1 + fabs (q) + fabs (data[i].desired_output.newq)) ||
              ( fprintf (stderr, "%g ?= %g\n", q, data[i].desired_output.newq),
                0 ));
      
      assert (fabs (logbeta - data[i].desired_output.logbeta) <=
              1e-4 * (1 + fabs (logbeta) + fabs (data[i].desired_output.logbeta)) ||
              ( fprintf (stderr, "(oops ... but this is ok) %g ?= %g\n", logbeta, data[i].desired_output.logbeta),
                1 ));

    }
}

int 
main (void)
{
  srand48 (69);

  test_mstep ();

  return 0;
}
